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ABSTRACT 

The orbital dynamics of most planetary satellites is governed by the quadrupole 
moment from the equatorial bulge of the host planet and the tidal field from 
the Sun. On the Laplace surface, the long-term orbital evolution driven by the 
combined effects of these forces is zero, so that orbits have a fixed orientation 
and shape. The "classical" Laplace surface is defined for circular orbits, and 
coincides with the planet's equator at small planetocentric distances and with 
its orbital plane at large distances. A dissipative circumplanetary disk should 
settle to this surface, and hence satellites formed from such a disk are likely to 
orbit in or near the classical Laplace surface. This paper studies the properties of 
Laplace surfaces. Our principal results are: (i) if the planetary obliquity exceeds 
68.875° there is a range of semimajor axes in which the classical Laplace surface 
is unstable; (ii) at some obliquities and planetocentric distances there is a distinct 
Laplace surface consisting of nested eccentric orbits, which bifurcates from the 
classical Laplace surface at the point where instability sets in; (ii) there is also 
a "polar" Laplace surface perpendicular to the line of nodes of the planetary 
equator on the planetary orbit; (iv) for circular orbits, the polar Laplace surface 
is stable at small planetocentric distances and unstable at large distances; (v) 
at the onset of instability this polar Laplace surface bifurcates into two polar 
Laplace surfaces composed of nested eccentric orbits. 

Subject headings: planets and satellites: formation - planets and satellites: gen- 
eral 
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Introduction 



In his study of Jupiter's satellites, iLaplacd (118051 ) recognized that the combined effects of 



the solar tide and the planet's oblateness induced a "proper" inclination in satellite orbits 
with respect to Jupiter's equator. He remarked that this proper inclination increases with 
the distance to the planet, and defines an orbital plane for circular orbits that lies between 
the orbital plane of the planet's motion round the Sun and its equator plane. This plane is 
called the Laplace plane. 

More generally, the Laplace plane is usually defined as the plane normal to the axis 
about which the pole of a satellite's orbit precesse£l The "Laplace surface" is the locus of 
all orbits that do not precess (i.e., the secular motion of the node and apse vanishes). 

In the most common situation, we consider circular satellite orbits around an oblate 
planet with non-zero obliquity, traveling around the Sun. The Laplace surface is then de- 
termined by the competition between the interior quadrupole potential from the equatorial 
bulge and the external quadrupole potential from the Sun. Close to the planet, the "classical" 
Laplace surface nearly coincides with the planetary equator, while at large distances it nearly 
coincides with the planetary orbital plane. The transition between these two orientations 
occurs near the "Laplace radius," defined below in equation 



The Laplace surface is important because it traces the shape expected for a thin gas 
disk or dissipative particulate ring surrounding the planet. Thus it is not surprising that 
many planetary satellites orbit close to the Laplace surface; those that do probably formed 
from a circumplanetary gas disk while those that do not were presumably either captured 
from heliocentric orbits or experienced unusual events in their past history. 

The purpose of this paper is to study the properties of the Laplace surface, including 
the stability of is generating orbits and its generalization to eccentric orbits. Although the 
Laplace surface has been known and studied for over two centuries, we believe that many of 
the results we present are new. 



1 Unfortunately, the term is sometimes also applied to the invariable plane, the plane perpendicular to 
the total angular momentum of an iV-body system. 
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2. Secular equations of motion 

2.1. The Hamiltonian 

The Kepler Hamiltonian that describes the motion of a test particle orbiting an isolated 
point mass M is 

h -w GM - GM - m 

Hk - ~ 2 v - — - - — , (1) 



r 



here r is the position vector measured from the center of the planet, v = r, r — |r|, and a 
is the semimajor axis of the test particle. The constants of motion are the Hamiltonian Hk 
(or semimajor axis a), and the angular momentum and eccentricity vectors 

1 r 

L = r x v, e = — v x (r X v) - -; (2) 

these are related to the eccentricity e and semimajor axis of the orbit by L 2 = GMa(l — e 2 ) 
and |e| = e. 

We now examine how the motion of the test particle is affected by additional forces 
from the equatorial bulge of the planet and the Sun. The quadrupole potential arising from 
an oblate planet is 

GMJ 2 R 2 GMJ 2 R 2 
*p00 = —^P^ose) = ^^[3(r • n p ) 2 - r 2 ], (3) 

where J 2 is the quadrupole gravitational harmonic, R p is the planetary radius, P 2 (x) = 
|(3x 2 — 1) is a Legendre polynomial, and 9 is the polar angle measured from the planet's 
spin axis, which is oriented along the unit vector n p . 

The quadrupole potential of the planet may be enhanced by inner satellites. If the planet 
hosts n satellites with masses rrii, i — 1, . . . ,n, on circular orbits in the equatorial plane of 
the planet with semimajor axes a,, then at distances r > the gravitational potential due 
to the satellites can be accounted for by augmenting J 2 to J 2 , where 

n 

J' 2 R\ = J 2 R\ + | a2 i m i/ M - ( 4 ) 

i=l 

Values of J 2 and J 2 for the giant planets and Pluto are given in Table [Q 

Since the solar tide is assumed to be weak, we may estimate its effects by averaging 
over the solar orbital period. We assume that the planetary orbit has semimajor axis a & and 
eccentricity e Q , and denote the normal to the orbit by the unit vector n Q . The obliquity of 
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the planet is then = cos 1 n ■ n p . Since r < a we need keep only the quadrupole term 
in the averaged solar potential, 

. GM Q 



4a|(l-e|)3/2 



[3(r-n e 



0; 



(5) 



The solar tide also causes the spin axis of the planet to precess, through its torque on the 
equatorial bulge. We neglect this effect to keep the analysis simple, since the precession rate 
of the planetary spin due to solar tides is n ormally much smaller than the precession rate of 
the satellite orbit (see iGoldreichl Il965l and lBoue fc Laskarll2006l for treatments that include 
precession of the planetary spin). 



The planet's "radius of influence" or "Hill radius" is 

M N 1/3 



3M 



(6) 



The Hill radius is roughly the point at which |$ | ~ GM/r; beyond the Hill radius the 
gravitational force experienced by a satellite is dominated by the solar tide rather than the 
force from the planet and most orbits are not bound to the pl anet. The Hill radius als o 
marks the location of the collinear Lagrange points of the planet ([Murray fc Dermottlll999l ). 
Our use of the orbit-averaged solar potential requires that the satellite radius r<rn. 



The total non-Keplerian potential due to the oblate planet, inner satellites, and Sun is 



Now average over the Keplerian orbit of the test particle, which has semimajor axis a, 
eccentricity e, and orientation specified by the unit vectors n along the angular momentum 
vector, u towards pericenter, and v = n x u. We have 



(r 2 ) 



<(r-u) 2 > 



a 2 (± + 2e 2 ) 



a 3 1 



e 2)3/2 • 



r ■ u 



r • v 



<(r-v) 2 > 

2 ' 



a 2 (| - |e 2 ) 



(7) 



2a 3 (1 



2)3/2' 



The averaged potential is $ = $ p + $ , where 



© 



GMJLRl 

' ' [l-3(n p -n) 2 ] 



4a 3 (1 



e 2)3/2 



3GM a 2 



4*1(1 



-0 



)3/2 



[(i + 2e 2 )(n -u) 2 + (I-ie 2 )(n G 



v) 2 - h 2 



The averaged Hamiltonian is then 



H 



H K + % + $ . 



|] • 0) 



(10) 
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By virtue of the orbit averaging, this Hamiltonian is independent of the mean anomaly, so 
its conjugate momentum (GMa) 1 ' 2 is a constant of motion, which in turn means that the 
semimajor axis may be treated clS db constant. 

Now set 



21/2 \GM M a 3 J' 2 R 2 p 

J = (1 " e) n ' e = 6U ' T= V^' e&= MaUl-elfl^ €p = ^- 

(11) 

The angular momentum L = (GMa) 1//2 j and, as suggested by the notation, e is the eccen- 
tricity vector (eq. [2]). If we define a dimensionless potential ^ p + \I/ = ($ p + $ )o/(GM), 
we have 

*p = 4(1 _ e ; 2)5/2 i 1 - e2 - 3 (j ■ n ^) 2 ] - 

*e = ^[5(e.n ) 2 -(j.n ) 2 -2e 2 ]; (12) 

an unimportant constant has been dropped. 

To repeat, the Hamiltonian + GM($ P + ^ Q )/a is based on the assumptions that 
(i) the precession rate of the planetary spin due to solar tides is negligible; (ii) the satellite 
is a massless test particle; (iii) the Sun is far enough from the planet that the solar tide can 
be approximated by a quadrupole; (iv) the satellite is far enough from the planet that the 
potential from the planet and the inner satellites can be approximated as a a monopole plus 
a quadrupole; (v) the perturbing forces due to Hf p + ^ Q are weak enough that the secular 
equations of motion can be used to describe its orbital evolution. 



2.2. Equations of motion 

Using equation (j2J) it is straightforward to show that the Poisson brackets of j and e are 

= v&m^ {ei) ej} = 7MR Eijkjk ' {jh ej} = vdm € * jlfih ' (13) 

where is the antisymmetric tensor. 

The time evolution of any variable / determined by the Hamiltonian H is given by 

!={/,*}• (14) 

In secular dynamics the semimajor axis is fixed and H can be considered to be a function 
only of the shape of the orbit, as expressed by j and e. Then from the chain rule 

^ = {/,j}Vj^+{/,e}V e ^ (15) 
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where Vj is the vector (d/dj\,d/dj2,d/dj$) with a similar definition for V e . Replacing / 
successively by ji and and using the relations f|T3|) . we find 



dt 
de 
~dt 



1 



VGMa 
1 

VGMa 



(ex V e # + jx ViH) 



(j x V e H + e x Vj#) . 



(16) 



Since e and j are constants of motion for the Kepler Hamiltonian Hk, the contribution of Hk 
to the right side of this equation must vanish, so we can replace H = Hk + ® by <fr. Equations 
ffTB"j) date back to iMilankovichl (119391 . eqs. 204 and 205; see also lBreiter fc Ratajczaki 12005 
and references therein). 

These equations admit three integrals of motion, $, je and j 2 +e 2 . Physically meaningful 
solutions are restricted to the four-dimensional manifold on which 



je = 0, 



•2 , 2 

J +e 



1. 



(17) 



Replacing t and $ by the dimensionless variables r and \I/ defined in equation ( flTl) we 



obtai 



Kcp 
rf 



J- = -j x Vj^ -ex V e f, 

dr 

d e _ T . _ T 

_ = -exVi$-ixV e i 
rfr j J e 



For the potential given by equations (fT2l we finally have 



4j 

de 
d^ 



3e j ■ n e 
4 

3e j • n 



J x n 



15e Q e ■ n P 



o 



e x n 







3fp j ■ n p 

2(1 _ e 2)5/2 



J x n. 



/j- 



4 

+ 



-e x n Q 



15e e-n . 3e p j ■ n p 
j x n + -TT- - — x n, 



2(1 - e 2)5/2 



3 3 1 - e 2 - 5(j ■ n p ) 

2 e — ~ e r. 



4 C f (1 _ e 2)7/2 



J x e. 



(19) 



To avoid distraction by trivial cases, we shall always assume that e p > (the planetary 
quadrupole is non-zero and positive, as expected for a planet that is oblate or has inner satel- 
lites), e > (solar perturbations are non-negligible), and n is neither parallel, antiparallel, 



2 There is a gauge freedom in the definition of ^(j,e) since the variables j and e are related by (fT7|) . ft is 
shown in Appendix [XI that the secular equations of motion (|18p are independent of gauge. 
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or perpendicular to n p (the planetary obliquity is not 0, ±|7r,7r). The plane defined by the 
planetary spin axis n p and the normal to the solar orbit n Q will be called the principal plane. 

We shall sometimes use cylindrical coordinates with the z-axis oriented along n p x n , so 
the plane z = coincides with the principal plane. The positive x-axis is chosen to coincide 
with n p . Then n Q lies in the z = plane, so it may be specified by its azimuthal angle <p & 
(the obliquity), which lies in the range (0,7r). 



In the limit e p — > equations (1191) provide a vector description of Kozai oscillations, the 
secular oscillations in the eccentricity and inclination of ( for example) a p lanet orbiting a 
member of a binary star ( Kozailll962l ; iHolman et al.lll997l ; iFord et al.ll200( 



Equations ( fl9l) are invariant under the transformations 



-e or rifT 



■0 



-n or n p 



-n„ or (j 



-r . 



(20) 



Invariance under n Q — > — n implies that we can restrict the range of the obliquity from 
(0,vr) to (0,|tt). 

Recall that equations (fl9l hold for arbitrary eccentricity, i.e., they do not represent an 
expansion that is valid only for e <C 1. 

We define the Laplace equilibria to be stationary solutions of equations (ITUl) . and the 
Laplace surface(s) to be the locus of all orbits that are Laplace equilibria. 



3. Circular Laplace equilibria 

The equations of motion (1191) can be solved explicitly in the case of circular orbits (e = 0). We 
banish this discussion to Appendix [B] and focus in this section on the equilibrium solutions, 
with e = and j =constant, which we call the circular Laplace equilibria. 

In circular Laplace equilibria the second of equations ffl9l is satisfied trivially. The first 
of equations f[T§j) yields 

e©(j ■ n )j x n + 2e p (j ■ n p )j x n p = 0. (21) 

Taking the scalar product of this equation with n Q , and again with n p , we conclude that 
either (i) j • n p = j ■ n = 0; (ii) j ■ (n p x n ) = 0. In the first case j is perpendicular 
to the principal plane, and we call this the "orthogonal" or "circular orthogonal" Laplace 
equilibrium; in the second case j lies in the principal plane and we call this the "coplanar" 
or "circular coplanar" Laplace equilibrium. 
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In the coplanar Laplace equilibrium, j lies in the z = plane, so it may be specified by 
its azimuthal angle 0. The equilibrium condition (121"]) becomes 

e sin 2(0 -0 ) + 2e p sin 20 = 0. (22) 

This equation has four solutions for in a 2ir interval; if is a solution then ± ^ir and 
+ 7r are also solutions. 

Equation (I2"2l may also be written as 

a 5 sin 2(0 -0 ) + 2rfsin20 = 0, (23) 

where the Laplace radius Tl is defined by 

rl = J' 2 Rla%{l-elf'^. (24) 



3.1. Stability 

The stability of the circular Laplace equilibria can be determined by writing j = jo + ji and 
expanding equations (|T9|) to first order in ji and e: 



dji 3e ji ■ n . 3e j o ■ n . 3e p ji ■ n p . 3e p j ■ n 



jo x n H — -ji x n H ^ — -j x n p H ^- — -ji x n 



dr 4 " " 4 

de 3e j o • n 15e e • n . 3e p j • n p 

— = exn j o xn + e x n p 

+ {f e - hpl 1 - 5 (jo • n P ) 2 ]} Jo x e. (25) 

Note that the two equations are decoupled: the linearized evolution of the orientation of the 
orbital plane, specified by j, is independent of the linearized evolution of the eccentricity and 
apse direction, specified by e. 

It is easy to show using the equilibrium equation (|2~T1) that a trivial solution of the first 
of these equations is ji = kjo where k -C 1 is a constant. In other words, the eigenvalue 
equation in A obtained by assuming ji oc exp(At) always has one zero eigenvalue. This is 
an unphysical solution of the linearized equations since the constant of motion j 2 + e 2 = 1 
requires that jo • ji = 0. 



3.1.1. The circular orthogonal Laplace equilibrium 

In the orthogonal equilibrium equations (|2"5"1) simplify to 

dji 3e ji • n . 3e p ji • n p . 

= 4 Jo x n © + 2 Jo x n p , 
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Fig. 1. — The regions in which stable, circular, coplanar Laplace equilibria exist are marked 
by heavy black stippling and the label "00". In the more lightly stippled regions, circular 
coplanar equilibria exist but are unstable, according to equations fl32l) and fl34l) . The first 
figure in each label is "0" or "1" according to whether the solution is stable or unstable to 
changes in the orbit plane orientation, described by j. The second figure is "0" or "1" ac- 
cording to whether the solution is stable or unstable to changes in the eccentricity, described 
by e. The vertical coordinate is the inclination of the orbit relative to the planetary equator, 
0. The horizontal coordinate is the obliquity of the planet relative to the ecliptic, <j>@. The 
results are unchanged if <fi — > <fi + tt or — > + 7r or (0 , 0) — > (n — (f> & , 7r — 0). 



-10- 



_e = _15eoe^i0 jo><n0 + ^ eQ _|g joXe (2g) 

These can be converted to eigenvalue equations by assuming ji,ei oc exp(Ai). To analyze 
the first equation, set u = ji ■ n , v = ji ■ n p . Taking the scalar product of dji/dt with n 
and n p we find 

du • dv 

— = §e p n • (j x n p )v ; — = f e n p • (j x n )w. (27) 

these can be combined to show that either (i) u — v — 0, A = 0; this is the unphysical 
solution noted above; or (ii) 

A 2 = -|e p e [j o -(n xn p )] 2 . (28) 

since A 2 < in this case, j x is oscillatory, so the circular orthogonal equilibrium is linearly 
stable to variations in the angular momentum vector j. The quantity in square brackets is 
just sin 00, where is the obliquity 

The second of equations (126]) can be analyzed similarly A trivial solution is e oc jo, 
A = 0; this solution is unphysical since the constant of motion j • e = requires jo ■ e = 0. 
The other eigenvalues are given by 

A 2 = ^(3e + e p )(2e -e p ). (29) 

Thus the circular orthogonal equilibrium is stable to variations in eccentricity e if and only 
if 

2e <e p or a < 2~ 1/5 r L . (30) 



3.1.2. The circular coplanar Laplace equilibrium 

In the circular coplanar equilibrium, the first of the linearized equations ( |25i) implies that 
ji oc exp(At) where 

= -j4 cos2< ^~ T6 e e cos2 (^ ~ 00 ) - ^e p e [cos20 + cos2(0-0 ) +2 cos 20©]. (31) 

Substituting the equilibrium condition (|22|) we have 

2 9e sin 2 cos 

A = ~ ifl ■ 2 7 — cos(0-0 ). (32) 

lb sin 0cos0 

The regions in which A 2 > (instability) are marked by red and blue stippling and the labels 
"10" and "11" in Figure [TJ All equilibria with < (ft < are stable to perturbations of this 
kind (variations in j but not e) and equilibria with < <f) < n are unstable. 
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The circular coplanar equilibria that are stable in this sense have — > as a — > and 
— ► as a — ► oo, so the corresponding Laplace surface coincides with the equator of the 
planet at small radii and with the planet's orbital plane at large radii. We call this the 
"classical" Laplace surface since this was the surface discovered by Laplace and the one that 
has been the focus of most work on this subject. 

The second of the linearized equations (|25|) implies that e oc exp(At) where 

9e 2 9f 2 

A 2 = v - [5 cos 4 - 2 cos 2 + 1] + — ^ [6 - 7 cos 2 (0 - Q )1 (33) 

16 16 

H — rr^[5 — 6 cos 2 — 6 cos 2 (0 — ) — 6 sin 20 sin 2(0 — ) — 9 cos 20 cos 2(0 — )1. 
64 

Substituting the equilibrium solution (|2"2"|) we have 
9e 2 r 

A 2 = ^ — <^ - 106 + 24 cos 20+ 146 cos 40- 100 cos(60 - 20 ) 

2048 sin 20 I 

-24cos(20 - 40 ) + 224cos(20- 20 ) - 54cos(40 - 40 ) - 8cos20 - llcos40 
-124cos(20 + 20 ) +25cos(80-40 ) + 8cos(40-20 )|. (34) 

The regions in which A 2 > (instability) are marked by green and blue stippling and the 
labels "01" and "11" in Figure U 

The inclination to the planetary equator is plotted as a function of the strength of 
the planetary quadrupole e p in Figure [2} for obliquities = 10°, 20°, . . . , 80°. Solutions for 
> 90° can be obtained by the transformation (0, ) — > (n — 0, 7r — ). 

Unstable equilibria are shown in Figure [2] by dotted or dashed lines. All equilibria with 
90° < < 180° are unstable. In addition, some of the classical equilibria (0 < 90°) are 
unstable to eccentricity growth. These are shown as stippled or shaded regions in Figures 
E] and |H Instability first appears at obliquity = 68.875° and is restricted to semimajor 
axes between about 0.9 and 1.25 times the Laplace radius r L (eq. [2] 



4. Eccentric Laplace equilibria 

We now look for stationary solutions to equations (fT9l) in which the eccentricity e = |e| is non- 
zero (recall that |j| = (1 — e 2 ) 1 / 2 ). It can be shown (see Appendix that all such solutions 
either have both j and e in the principal plane defined by n and n p (the "coplanar-coplanar" 
or "eccentric coplanar-coplanar" solution), or one of j and e in the principal plane and the 
other orthogonal to it (the "coplanar-orthogonal" or "orthogonal-coplanar" equilibrium if j 
or e, respectively, lies in the principal plane). 
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0.2 0.4 0.6 0.8 1 

e o/( e p + e ©) 

Fig. 2. — Circular coplanar Laplace equilibria. The vertical axis is the angle <fi between the 
pole of the planet and the orbital pole of the satellite (the inclination of the orbit relative 
to the planet's equator). Solutions are shown for eight values of the planetary obliquity, 
00 = 10°, . . . , 80°. There are two equilibrium curves for each obliquity, one with < <p < 90° 
and the other with 90° < <fi < 180°. Solid lines denote stable equilibria, while dashed and 
dotted lines denote equilibria with one or two unstable roots respectively. All equilibria 
with <p > 90° are unstable. The classical equilibria are those with < 90°. The horizontal 
axis represents the relative strength of perturbations from the solar tide and the planetary 
quadrupole (eq. ITT1) ; the planetary quadrupole dominates on the left side of the figure and 
the solar quadrupole on the right. 
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Fig. 3. — Regions in which the classical Laplace equilibria are unstable are stippled. The 
top panel shows the unstable range of obliquities as a function of the ratio of the solar 
and planetary perturbation strengths (eq. [11]) and the bottom panel shows the unstable 
obliquities as a function of semimajor axis (in units of the Laplace radius r^, eq. |2^|) . All 
instabilities are in the eccentricity vector; the angular-momentum vector is stable. 
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Fig. 4. — Regions in which the classical Laplace equilibria are unstable are shaded. 
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4.1. Eccentric coplanar-coplanar Laplace equilibrium 

In this case j and e lie in the principal plane, and equations f|T9|) with dj/dr = de/dr = 
yield 

2e p sin20 = e (l - e 2 ) 3/2 (l + 4e 2 ) sin 2(0 - 0) 
e p [l-3cos 2 0] = e (l-e 2 ) 5/2 [l-4sin 2 (0 -0)]. (35) 

These can be used to solve for e 2 and 0, given and e p /e . The solutions are physical only 
if < e 2 < 1. 

The regions in -0 space in which eccentric coplanar-coplanar Laplace equilibria exist 
are shown in Figure by stipples. These regions are bounded, in part, by the lines = , 
= + and = ± cos _1 (l/v/3) = 54.7°, 125.3°. Heavy stippling marks regions in 
which the equilibria are stable. 

The stable regions are small for > 54.7° so we will focus on the region < 54.7°. The 
eccentric coplanar-coplanar equilibria in this region are closely related to the circular coplanar 
equilibria discussed in £J3J We found there that the circular coplanar Laplace equilibria were 
unstable for a range of e p /e when the obliquity > 68.875°. Figure [6] shows the range in 
which the circular coplanar equilibria are stable for = 70°, 71°, . . . , 89° as solid horizontal 
lines, with gaps marking the unstable range. Superimposed on these lines are the eccentric 
coplanar-coplanar equilibria, marked by heavy blue or light red curves depending on whether 
they are stable or unstable. The height of these curves is proportional to the eccentricity of 
the equilibrium. The figure shows that the eccentric coplanar-coplanar equilibria bifurcate 
from the circular coplanar equilibrium at the point where the circular coplanar equilibrium 
becomes unstable. The structure of the secular Hamiltonian at the bifurcation is that of 
the standard resonant Hami ltonian at j + 2 : j orbital resonances (k = 2 in the notation of 



Borderies k, Goldrei chill 9841 ) . 



For 68.875° < < 71.072° the eccentric coplanar-coplanar solution exists in a limited 
range of e p /e & and is stable throughout this region. For > 71.072° the eccentric coplanar- 
coplanar solution is unstable for some part of the range of e p /e in which it exists. Let us 
imagine a satellite in the circular coplanar Laplace equilibrium at (say) e /(e p + e ) = 0.8. If 
e /(e p +e ) slowly decreases (for example, because the satellite is slowly migrating toward the 
planet) then we expect that (i) for < 68.875° the satellite will always remain in a circular 
coplanar Laplace equilibrium; (ii) for 68.875° < < 71.072° the satellite will transfer onto 
the eccentric coplanar-coplanar equilibrium, its eccentricity will then grow as its spirals in, 
reach a maximum depending on , then shrink back to zero, at which point it will rejoin 
the sequence of circular coplanar Laplace equilibria and spiral into the planet on a circular 
orbit; (iii) for > 71.072° the satellite will transfer onto the eccentric coplanar-coplanar 
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equilibrium, and its eccentricity will grow until the equilibrium orbit becomes unstable, 
at which point it presumably undergoes large oscillations in eccentricity and inclination. 
Numerical simulations of this evolution are described in §3 we shall find that (i) is correct, 
but that (ii) and (iii) need to be qualified in that orbits track the sequence of eccentric 
equilibria less well than they imply. 

An interesting but difficult question, which we do not attempt to answer here, is the 
behavior of a dissipative gas or particulate disk in the region where the classical Laplace 
surface is unstable. There are at least three alternatives: (i) The dissipative forces suppress 
the secular instabilities in the circular coplanar Laplace equilibria, thereby allowing the disk 
to occupy the classical Laplace surface, (ii) The disk occupies the surface defined by the 
eccentric coplanar-coplanar equilibria. A possible problem is the large eccentricity gradients 
in this sequence of equilibria: a flat, cold disk composed of aligned eccentric orbits cannot 
exist if \a(de/da)\ > 1, because the orbits cross. This condition does not apply directly to the 
Laplace surface because it is not flat; nevertheless, the large eccentricity gradients may lead 
to strong shears that destabilize the disk, (iii) A gap opens in the disk where the classical 
Laplace surface is unstable. 



4.2. Eccentric coplanar-orthogonal Laplace equilibrium 

In this case j lies in the principal plane and e is normal to it, so 

n Q • (j x n p ) = 0, e • n = 0, e • n p = 0. 

The first of equations f[T§|) with dj/dr = yields 

e (l - e 2 ) 5/2 sin 2(0 - ) + 2e p sin 20 = 0. 



(36) 



We introduce a vector w defined bywxj = eorw = jxe/(l — e 2 );w lies in the principal 
plane and w- j = w - e = 0. Substituting for e in the second of equations (fT^j) with de/dr = 
and e • n = yields 







e j-n [j(w-n )-w(j-n J 



2e P j • n p 

(1 - e 2 ) 5 /2 



[j(w • n p ) - w(j • n p )] 



+ 



2e m 



5(j • % 



fw. 



(37) 



-® C P (1 _ e 2)7/2 

The vectors j and w are linearly independent, so the coefficient of each must be zero. The 
coefficient of j vanishes if equation fl36|) is satisfied. The condition that the coefficient of w 
vanishes is 



e [2 — cos 2 



+ 



(1 - e 2 ) 5 /2 



(3 cos 2 



0. 



(3? 
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Fig. 6. — The relation of the circular coplanar Laplace equilibria to the eccentric coplanar- 
coplanar equilibria. The gaps in the solid horizontal lines represent regions in which the 
circular coplanar equilibria are unstable. The curved lines represent the eccentric coplanar- 
coplanar Laplace equilibria with inclination < 90°. The y-coordinate of these lines is 
00 + 20e where e is the eccentricity of the solution; the curves are heavy blue or light red 
according to whether the solution is stable or unstable. The figure shows that the eccentric 
coplanar-coplanar equilibria bifurcate from the circular coplanar equilibria. 
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Equations fl36l) and fl38l) can be combined to eliminate e, e , and e p : 

2[2-cos 2 (0 -0)]sin20+ (3cos 2 0- l)sin2(0 -0) = 0. (39) 

This result determines the inclination between the satellite orbit and the planetary spin in 
terms of the obliquity , independent of e p , e Q , or e. Numerical solution of this equation for 
< 00 < \k shows that there are two solutions for each value of Q , but the upper solution 
is unphysical, because these values of and yield e (l — e 2 ) 5 ^ 2 /e p < in equation ( l36l) . 
Thus there is a unique inclination for each value of the obliquity , as shown in the top 
panel of Figure H 

The bottom panel of Figure [7] shows the eccentricity of the coplanar-orthogonal Laplace 
equilibria. All of these equilibria are unstable. 



4.3. Eccentric orthogonal-coplanar Laplace equilibrium 

In this case e lies in the principal plane and j is normal to it, so 

j • n p = j ■ n = 0, e ■ (n p x n ) = 0. (40) 

The solution of equations ffl9|) with dj/dr = de/dr = requires that 

e-n = O, 2e (l-e 2 ) 5 / 2 = e p . (41) 

Solutions exist whenever 2e > e p . 

Linear stability analysis shows that small perturbations in j or e grow as exp(At) where 

q 8/5 2/5 99c, 

A 2 = ^Fsin 2 or A 2 = [1 _ ( Cp /2e ) 2/5 ] ; (42) 

since A 2 < the eccentric orthogonal-coplanar solutions are stable. 

Comparison with the results of §3.1.11 shows that these (eccentric) equilibria bifurcate 
from the (circular) orthogonal equilibrium sequence at the semimajor axis a = 2~ l ^ri where 
the latter becomes unstable. 



5. Orbital evolution 



Consider the evolution of a satellite on a circular orbit that is slowly decaying, so the satellite 
is migrating toward the planet. Let us assume that the initial semimajor axis is much larger 
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orthogonal Laplace equilibrium, as obtained by solving equation (|39l) . (Bottom) Eccentricity 
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Fig. 8. — Numerical integrations of the equations of motion (fl9l) . The ^/-coordinate is Q +2Oe 
where e is the eccentricity of the solution. Each orbit is begun at semimajor axis a = 1.3rt 
and the semimajor axis then shrinks according to a/ri = 1.3 — 0.00 2r where r is defined in 
equation fill I) . The satellite is initially in the circular Laplace equilibrium except for a seed 
eccentricity of 0.0001. The x-coordinate is e /(e p + e ) = a 5 /(a 5 + r£). Four integrations 
are shown, for obliquities = 70°, 75°, 80°, 85°. 
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than the Laplace radius ri (eq. [241 and that the orbital plane coincides with the planetary 
orbit (i.e., the satellite is in the classical Laplace surface). As outlined at the end of §4.11 if the 
obliquity <p Q < 68.875° we expect that the satellite will remain in a circular coplanar Laplace 
equilibrium as it spirals in. For 68.875° < < 71.072° we expect that at semimajor axis 
a ~ 1.2rx (Fig. [3]) the satellite will transfer onto the eccentric coplanar-coplanar equilibrium; 
its eccentricity will then grow as it migrates inward, reach a maximum depending on Q , 
then shrink back to zero, at which point (a ~ 0.9r^) it will rejoin the classical Laplace 
surface and spiral into the planet on a circular orbit. Finally, for Q > 71.072° the satellite 
will transfer onto the eccentric coplanar-coplanar equilibrium, and its eccentricity will grow 
until its orbit becomes unstable and begins to execute large oscillations in eccentricity and 
inclination 

Numerical integrations of the equations of motion ffH?|) are shown for a migrating satellite 
in Figure [H As expected, when the initial obliquity 0© is 75°, 80°, or 85° large eccentricity 
oscillations develop in the region where the circular coplanar Laplace equilibrium is unstable 
(compare Fig. [6]). Smaller oscillations develop in the case 0© = 70°; these are unexpected 
since the eccentric coplanar-coplanar sequence is stable at this obliquity, so we would expect 
smooth growth and decline in the eccentricity as the semimajor axis declines, as in Figure El 
Further integrations show that the behavior of migrating satellites at this obliquity depends 
on the migration time and other parameters; for example, the orbit shown in Figure [9] follows 
the eccentric coplanar-coplanar sequence for a while and then rather suddenly jumps to an 
orbit with a chaotic appearance that reaches eccentricities as high as 0.45. These results 
suggest the presence of large chaotic regions and perhaps higher order resonances in the 
phase space, but we have not yet explored these features. 



6. Applications 

Direct applications of these results to the solar system are rather few. Only Uranus and 
Pluto have obliquities that exceed the critical value 0© = 68.875° at which the classical 
Laplace surface becomes unstable, and these do not have satellites close to the unstable 
range of semimajor axis. 

The circular orthogonal Lapl ace equilibria have been suggested as possible sites for 



polar" rings around Neptune (see iDobrovolskis et al.lll989l and references therein). In this 



case the solar quadrupole potential is much weaker than the potential from Neptune's satellite 
Triton. However, the effects of Triton can be modeled approximately using the formalism 
we have derived, simply replacing the Sun by Triton. We find from equation (124]) that the 
Neptune- Triton Laplace radius is rt = 2.15 x 10 10 cm = 8.54R p . Equation (1501 then implies 
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Fig. 9. — Numerical integration of the equations of motion (fl9l) as a satellite migrates 
inward. The obliquity is = 70° and the satellite semimajor axis shrinks according to 
a/r^ = 1.112 — 0.000002r where r is defined in equation ffTTT) . The satellite is initially in 
the circular Laplace equilibrium except for a seed eccentricity of 10~ 6 . The x-coordinate 
is e Q /(e p + e ) = a 5 /(a 5 + r\). The satellite initially lies on the classical Laplace surface; 
at e Q /(e p + e ) = 0.621 it transitions to the eccentric coplanar-coplanar sequence; and at 
e o/( e P + e o) — 0.593 it transitions to an orbit that exhibits large and irregular eccentricity 
oscillations. 
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that circular polar rings would be unstable at semimajor axes > 7AR P , and beyond this 
radius polar rings must be eccentric. This analysis is only approximate, since (i) the Laplace 
radius is 61% of Triton's semimajor axis, a? = 3.548 x 10 10 cm, so the approximation 
of Triton's potential by its quadrupole component is poor; (ii) Triton's orbital axis and 
Neptune's spin axis both precess around their mutual invariable plane, whereas we have 
assumed that they are fixed; this is probably only a minor correction, since the precession 
period of 680 y is much longer than the growth time of the instability, (a 3 / GMp) 1 ^ 2 / A = 17.3 y 
at a = r L ; (iii) our calculation neglects collective effects in the ring, such as interparticle 
collisions, which might suppress the instability. 

Our original motivation for examining this problem was the orbit of Saturn's satellite 
Iapetus, which has an eccentricity of 0.028 and an inclination of 7.5° to the Laplace surface. 
If it formed from a circumplanetary disk, one might expect Iapetus to have zero eccentricity 
and inclination relative to this surface. Thus it appears that some proce ss has pumped up 
Iapetus's inclination while leaving its eccentricity near zero. IWardl (119811 ) pointed out that 
the shape of the classical Laplace surface is affected by the mass in the circumplanetary disk, 
and suggested that the current orbit of Iapetus reflects its shape before the disk dispersed. 
However, this scenario requires the dispersal of the disk in < 10 3 yr — if the dispersal were 
slower, the inclination relative to the Laplace surface would be an adiabatic invariant and 
thus would remain near zero. The semimajor axis of Iapetus a = 59i? p is not far from the 
Laplace radius tl = 48i? p (Table 1) so it is natural to wonder whether instabilities have 
excited the inclination. However, (i) at Saturn's obliquity of 26.7°, all circular orbits in the 
Laplace surface are stable; (ii) the instability we have found in the Laplace surface at high 
obliquity tends to excite eccentricity, not inclination. 

The rich dynamics described in this paper is likely to play a larger role in scenarios in 
which the obliquities of the giant plan ets have changed substantially since the formation of 
the solar system (e.g.. lTremaindll99ll ). 



A further application is to extrasolar planetary systems. It is well-known that Kozai 
oscillati ons can excite the eccentricity of a single planet orb iting one member of a binary star 
system (IKozail Il962l ; iHolman et al.l 119971 ; iFord et al.l 120001 ) . Now con sider a star that ho sts 
two planets, one a hot Jupiter, and belongs to a binary system (see iTakeda et al.l 120081 for 
numerical simulations of such systems). The hot Jupiter augments the (otherwise negligible) 
quadrupole moment of the host star, just as inner satellites augment the quadrupole moment 
of a planet (eq. @J. The resulting Laplace radius (I24"l) is 

1/5 



r L = 2.7 AU 



mj 



0.001M* 



2/5 



0.1 AU/ V300AU/ 



3/5 



1 - elf' 10 , 



(43) 



where mj and ctj are the mass and semimajor axis of the hot Jupiter, and M*, a*, and 
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are the mass, semi-major axis, and eccentricity of the binary companion. If, as we expect, 
the orbital axes of the hot Jupiter and the distant companion star are uncorrelated, 36% of 
binary stars will have obliquities that exceed the critical value 68.875° at which instabilities 
set in at some semimajor axis; and if planetary migration is common we expect that many 
planets now in the habitable zone will have passed through the unstable region and thereby 
acquired substantial eccentricities (cf. Fig. |SJ). 

7. Summary 

We have examined the orbital dynamics of planetary satellites under the combined influence 
of the quadrupole moment from the planet's equatorial bulge and the tidal field from the 
Sun. The Laplace equilibria are orbits in which the secular evolution due to these forces is 
zero. They represent the orbits in which planetary satellites formed from a circumplanetary 
gas disk should be found. 

Laplace equilibria exist for both circular and eccentric orbits. At a given semimajor 
axis, the orbit normals to the circular Laplace equilibria lie along three orthogonal directions, 
two of them in the principal plane defined by the planet's spin axis and the normal to its 
orbit around the Sun (the "circular coplanar" equilibria). One of the two circular coplanar 
equilibria is always unstable. The warped surface swept out by the other circular coplanar 
equilibrium orbits as the semimajor axis varies is called the classical Laplace surface (eq. 123]) . 
The classical Laplace surface coincides with the equator of the planet at small distances and 
with the orbital plane of the planet at large distances. Orbits in the classical Laplace surface 
are stable if the planetary obliquity Q < 68.875°, while in the range 68.875° < < 90° 
there is a range of semimajor axes near the Laplace radius ri (eq. |2"3|) in which the surface 
is unstable. The third Laplace equilibrium for circular orbits (the "orthogonal" or "polar" 
equilibrium) corresponds to orbits that cross over the planet's pole, and these are stable if 
and only if the semimajor axis is less than 2~ 1 ^ 5 ri (eq. [301) . 

At some obliquities and semimajor axes eccentric Laplace equilibria also exist. In the 
"coplanar-coplanar" equilibria both the eccentricity and angular-momentum vectors lie in 
the principal plane (Figures [5] and [6]). The "orthogonal-coplanar" equilibria are stable, 
eccentric polar orbits that bifurcate from the circular orthogonal equilibria at the semimajor 
axis 2~ 1 / 5 r/ j where these become unstable. Other eccentric equilibria exist but are unstable. 

The use of the secular equations of motion ffl9l should be legitimate so long as the 
precession times for the eccentricity and angular-momentum vector are much longer than 
the orbital period of the satellite. A sufficient condition for this is that the satellite semi- 
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major axis is much larger than the planetary radius and much smaller than the planet's Hill 
radius (cf. Tabled]). A possible limitation of the secular equations is that they neglect the 
evection resonance, where the apsidal precession rate equals the mean motion of the planet 
around the Sun. F o r the outer planets, the evection resonance typically occurs at ~ 0.2rz,- 



Touma fc Wisdom! (119981 ) have stressed the role of the evection resonance in the history of 



the lunar orbit. 

Several unanswered questions remain: (i) What is the structure of the four- dimensional 
phase space of the equations of motion (JT9"|) . and what fraction of the orbits are chaotic? 
These issues could be explored through surfaces of section, although the exploration would 
be laborious since there is a separate surface of section for each value of the Hamiltonian, 
semimajor axis, and planetary obliquity, (ii) What is the nature of the evolution of satellites 
in the classical Laplace surface that migrate into an unstable region (Fig. |Sj)? (iii) What 
are the shape and properties of dissipative disks in the range of semimajor axes where the 
classical Laplace surface is unstable? Does the dissipation stabilize the classical Laplace 
surface? Does the disk follow the eccentric coplanar-coplanar Laplace surface? Or is there 
no steady-state disk? Hints that eccentric structures can occur in dissipative disks include 
the eccentric structures seen in the planetary ring systems of Uranus and Saturn, debris disks 
around young stars, and galactic nuclei such as M31. (iv) What is the analog of the Laplace 
surface in accretion disks around black holes, where Lense-Thirring precession rather than 
a quadrupole potential is the dominant non-Keplerian perturbation from the central body? 

We thank Yuri Levin for thoughtful comments. This research was supported in part by 
NASA grants NNX08AH83G and NNX08AH24G, and by NSF grant AST-0807432. 



A. Gauge dependence in the averaged potential 

The orbit-averaged potential e) is a function of the six scalar components of j and e but 
is only physically meaningful on the four-dimensional manifold defined by the constraints 
( JT7|) . Thus \I/ may be replaced by ^ + F, where the gauge function F(j,e) is arbitrary 
except that F = const on the manifold (fTTI) . We now show that F(j, e) has no effect on the 
equations of motion (ITS]) , just as adding a constant to a Newtonian potential has no effect 
on the equations of motion x = — V0. 

We work in the coordinates with axes parallel to the unit vectors (u, v, n) defined just 
before equation ((7j). Substituting F for \l> in the equations of motion (ITBl and observing 
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that j = jn, e = eu, we have 



= _j v ( VjF • u) + ju(VjF ■ v) - en(V e F ■ v) + ev(V e F • n) 

dr 
de 

— — = -jv(V B F ■ u) + ju(V e -F • v) - en(VjF • v) + ev(V;F • n). 

dr 



(Al) 



The gauge condition that F = const on the manifold (|17|) is 

= dF = VjF ■ dj + V e F ■ de when e • dj + j ■ de = and j • dj + e ■ de = 0. (A2) 
We write <ij and de in terms of the unit vectors (u, v, n) as 

d] = a u u + a v \ + a n ii, de = b u u + b v v + b n n; (A3) 
then (|A2p requires that 

= a u { VjF ■ u) + a„( VjF ■ v) + a n { VjF • n) + b u { V e F ■ u) + b v { V e F • v) + 6 n ( V e F ■ n) (A4) 
when 

ea u + jb n = 0, ja n + eb u = 0. (A5) 

The six variables a u ,a v ,a n ,b u ,b v ,b n can be chosen arbitrarily so long as they satisfy the 
two constraints (IA5j) . Thus (IA4I) can only be true if the gauge function F satisfies four 
conditions: 



V,F • v = 0, V e F ■ v = 0, j V,F ■ u = eV e F ■ n, eVjF ■ n = j V e F ■ u. 



(A6) 



With these conditions, equations flAip yield dj/dr = 0, de/dr = 0, so the gauge function 
has no effect on the equations of motion. 



B. Behavior of circular orbits 

A circular orbit remains circular under the perturbations in question; that is, there is a set 
of solutions of equations (fTUj) in which e = for all time. The dynamics of the circular orbit 
is governed by the three-dimensional system of equations 

# = 3£Qj nQ j x n + ^^j x ry (Bl) 
dr 4 2 K J 

This system has two integrals of motion, | j | (equal to unity on the manifold corresponding 
to physical solutions) and the associated Hamiltonian (eq. [lOj for e = 0). The system is 
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integrable, as we now demonstrate by recasting the problem in the form of free rigid-body 
dynamics with the appropriate inertia tensor. 

First, note that (j -n )j x n = j x (j ■ n )n can be rewritten in the form j x T • j 
where 

(n\ n x n 2 nin 3 \ 
n 2 n x n\ n 2 n 3 (B2) 
n 3 rii n 3 n 2 nj / q 

is a symmetric tensor built out of the components of n = (711,712,713)0. The same can 
be done for the n p term leading to an equivalent symmetric tensor T p , built out of the 
components of n p . We can thus rewrite equation (1B1I) in the suggestive form 

A = j x T • j where T = aT & + (3T P , (B3) 

with a = |e and (3 = ~e p . Following this rearrangement, the program is straightforward: 
we find the matrix that diagonalizes the symmetric tensor T, identifying principal directions 
(actually the directions of the circular Laplace equilibria discussed in §3]), and the principal 
values, which will decide the stability of equilibria and the global phase-space topology 
around them. To fix things, it is simplest to take n p = (1,0,0), and n p = (cos0 , sin0 , 0) 
(obliquity ), leading to: 

a cos 4> Q sin (f) Q 

asin 2 I . (B4) 




The eigenvalues of T are 



h = 

a + P- V( a + P) 2 ~ sin2 



2 



+ P + y/(a + P) 2 - 4a/3 sin 2 



a 



(B5) 

with t\ = < t 2 < t 3 (i 2 = for (f> Q = 0, an ignorable case). We can then solve for 
the eigenvectors t^, t^ 2 \ of T, the principal axis directions. Since T is symmetric its 
eigenvectors are orthogonal. Moreover the eigenvectors can be chosen to be orthonormal and 
to form a right-handed triad (t^ = x t^). In particular, = (0,0, ±1) (assuming, 
as usual, that e p , e , and sin0 are non-zero), that is, t^ 1 ) is perpendicular to the principal 
plane containing n p and n , while t'- 2 -' and lie in the principal plane. 

Let R be the orthogonal transformation associated with these vectors, the transfor- 
mation that diagonalizes T. The columns of R are the orthonormal eigenvectors of T, so 
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R i3 = t[ j) . We now rotate to the principal- axis coordinate system, by setting j = R • J. 
Equation (1B3j) turns into 

dJi 

~~dt 

dJ 2 

~dT 

dJ 3 
dt 



(*8 






(*1 




— —t-zJiJz 


(*2 


-h)J x J 2 


= t 2 J\J 2 ^ 



(B6) 



exactly what one would obtain for the asymmetr ical top in Euler's formu lation, with 1/t 



standing for Jj, the principal moments of inertia (lLandau &: Lifshitzlll976l ). As in the case 
of the rigid body, vectors that lie along the principal axes (parallel and anti-parallel) are 
stationary states of the dynamics: J = (±1, 0, 0), (0, ±1,0), and (0, 0, ±1). 

The equilibria along t^ 1 ' are linearly stable, with linear oscillation frequency fl% = £3 
{t 2 and t 3 are both positive); equilibria along t^ 3 ^ are also stable, with linear oscillation 
frequency f2 2 = a/ (£3 — t 2 ) £3 (i 3 > t 2 > 0); equilibria along t^ 2 ^ are linearly unstable with 
growth rate ^2(^3 — t 2 ). Thus and t^ 2 ^ describe the directions of the angular-momentum 
vectors for the linearly stable and unstable circular coplanar Laplace equilibria, while 
describes the circular orthogonal Laplace equilibrium, which is stable (recall that we consider 
only circular orbits in this Appendix, so the stability properties described here do not include 
the eccentricity instabilities discussed in §3.11) . 

We now examine the nonlinear dynamics. The Hamiltonian is the restriction of ^/ p + ^q 
(eq. [T2l) to circular orbits, e = 0. In the present notation this can be written 

3 

H c = -|j T • T • j = -I = + ^ 3 2 )> (B7) 

i=i 

where the superscript "T" denotes transpose. If it is not clear already, note that H c generates 
Lie-Poisson dynamics of J according to: 

7 = -Jx V 3 H C . (B8) 

GST 

There are two integrals of motion, H c itself and the mag nitude J 2 = J\ + Jj + Jf = 1. 
Thus trajectories lie on the intersection of the angular momentum sphere and the elliptical 
cylinder representing the energy surface (as opposed to the ellipsoid with axes that is 
involved in the solution of the asymmetrical top) . The largest allowable energy is H c max = 
and the smallest is H Ctmin = —^3. Starting at the minimum energy, the cylinder intersects 
the sphere at two points, corresponding to equilibria along t^. Increasing the energy, the 
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cylinder intersects the sphere in closed curves around these equilibria. This remains so, till 
we reach the critical energy H c = — \t 2) at which point the cylinder is tangent to the sphere 
at ±t^ 2 ) as it intersects it along the limiting curves, separatrices, that separate the stable 
librations around t^ 3 ^ from the stable librations around t^ 1 ). Increasing further, the cylinder 
now intersects the sphere at two librating curves around ±t^\ until it reduces to a needle 
piercing the sphere at ±t^ when H c = 0. 

Orbit shapes are easy to obtain in projection. For trajectories librating around the 
equilibrium, take the ratio of the J\ and J 2 equations and integrate (or eliminate Jf between 
the energy and angular-momentum integrals) to get 

{h-t 2 )Jl +hJl = const = t 3 + 2H c , -t 3 <2H c <-t 2 . (B9) 

At the limiting value, H c = — \t 2) this defines the separatrices in projection. Similarly we 
can obtain trajectories librating around the equilibrium. Of course the equations of 
motion (|B6D can be solve d explicitly, following the methods used for the asymmetrical top 



( ILandau &: Lifshitzlll9761 ). 



Incidentally, if the system is subject to a process that dissipates energy while conserving 
angular momentum, it will decay towards the minimum energy state, which lies along ±t®, 
the direction normal to the circular coplanar Laplace equilibrium. This is often interpreted 
to mean that the classical Laplace surface, at a minimum of H c , is secularly stable while the 
circular orthogonal Laplace equilibrium at a maximum of H c is secularly unstable. However, 
most dissipative processes also affect the eccentricity and semimajor axes of the orbits, so the 
evolution of the orientation of th e angular-momentum vect or due to dissipation cannot be 



treated in isolation. In particular, iDobrovolskis et al.l (Il989l ) argue in the context of possible 



polar rings around Neptune that the circular orthogonal Laplace equilibria are secularly 
stable. 



C. Classification of eccentric Laplace equilibria 

We investigate the properties of Laplace equilibria (stationary solutions of eqs. [T9|) with non- 
zero eccentricity. As usual we assume that e p > and that n is neither parallel, antiparallel, 
nor perpendicular to n p . Recall that |j| = (1 — e 2 ) 1 / 2 . 

Take the scalar product of the first of equations (TT9l) with n . The first two terms on 
the right side vanish, and since dj/dr — in a stationary solution and e p ^ we must have 
either 

(a) n • (j x n p ) = 0, or (b) j ■ n p = 0. (CI) 
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First consider ). Take the scalar product of the second of equations ffT9l) with 

e. The first, third, and fourth terms on the right side vanish, and e 7^ 0, so we must have 
either 

(aa) n ■ (j x n p ) = and e ■ (j x n ) = 0, or 

(ab) n • (j x rip) = and e ■ n = 0. (C2) 

Case (aa) requires that j, e, n , and n p are all coplanar, that is, that j and e lie in the 
principal plane defined by n and n p . We call this the "coplanar-coplanar" solution since 
j and e lie in the principal plane. Case (ab) requires that j, n , and n p are coplanar and 
that e is perpendicular to n ; moreover by definition e is perpendicular to j. Hence j lies 
in the principal plane and e is normal to this plane; we call this the "coplanar-orthogonal" 
solution. 

Now consider case (b). Take the scalar product of the first of equations ffl9l) with e. 
The second and third terms on the right side vanish, so either 

(ba) j ■ rip = and j • n = 0, or 

(bb) j • np = and e • (j x n ) = 0. (C3) 
In case (ba), the first of equations ffTU]) implies that (e • n )e x n = 0. Hence either 

(baa) j • n p = and j • n = and e • n = 0, or 

(bab) j • n p = and j • n = and e x n = 0. (C4) 

In case (baa), j is perpendicular to the principal plane, and e lies in the principal plane at 
right angles to n ; we call this the "orthogonal-coplanar" solution. In case (bab) e is parallel 
or antiparallel to n so we can write e = ±en and the second of equations (TT^|) reduces to 

3e + (1 _ £ e P 2)5/2 = 0, (C5) 

which has no solution since e , e p > 0. 

In case (bb), e, j, and n lie in the same plane, and the first of equations (fT9l) reduces 

to 

(j • n )j x n = 5(e • n )e x n . (C6) 

Let e and j define the positive x- and y-axes of a Cartesian coordinate system and in this 
system let n = (cos ip, sin ip). Then equation flC6j) becomes (l + 4e 2 ) sin tp cos ip = 0. Hence 
either sinip = or cosip = 0. If simp = 0, e is parallel or antiparallel to n , and we return 
to case (bab), which has no solution. If cosip = 0, then j is parallel or antiparallel to n , so 
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n ■ n p = 0, a special case (planetary obliquity of 90°) that we have already excluded from 
consideration. 

Thus the only eccentric Laplace equilibria are case (aa), in which both j and e lie in the 
principal plane formed by n p and n (coplanar-coplanar equilibrium); case (ab), in which j 
lies in the principal plane and e is normal to this plane (coplanar-orthogonal equilibrium); 
or case (baa) , in which e lies in the principal plane and j is normal to the plane (orthogonal- 
coplanar equilibrium). 
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Table 1. Properties of the outer planets 



planet 


a Q (AU) 


Rp (km) 


J 2 


A 


obliquity 


ri/Rp 


r H /R P 


Jupiter 


5.2029 


71492 


0.014696 


0.045020 


3.1° 


35.36 


743.3 


Saturn 


9.5367 


60330 


0.016291 


0.070561 


26.7° 


48.40 


1080.1 


Uranus 


19.189 


26200 


0.003343 


0.018699 


97.9° 


63.96 


2675.1 


Neptune 


30.070 


25225 


0.00341 


0.024069 


29.6° 


93.20 


4600.8 


Pluto 


39.482 


1151 




14.296 


112.5° 


419.6 


6935.8 



Note. - - The planet's semimajor axis and radius are a Q and R p . Obliquity is the angle 
between the planet's spin and orbital axes, which in the notation of this paper is given by 
00 = cos -1 rip • n Q . The Laplace radius ri and Hill radius r# are defined by equations (1241) 
and ®. 



